Enhanced disk actuator modeling

ABSTRACT

Methods, systems, and apparatus, including computer programs encoded on computer storage media, for modeling turbine parameters. One of the methods includes obtaining, along multiple points of a blade of a turbine from a minimum radius rmin of the blade to a maximum radius rmax of the blade, lift coefficients C yi  and drag coefficients C xi . At the multiple points of the blade from rmin to rmax, corresponding components of an upstream fluid flow velocity vector u h,Ri  and u φ,Ri  and components of a downstream fluid flow velocity u h,Li  and u φ,Li  are obtained. Averaged directions β i  of the upstream and downstream fluid flow velocity vectors are computed using the components of the upstream fluid flow velocity vector u h,Ri  and u φ,Ri  and the components of the downstream fluid flow velocity u h,Li  and u φ,Li . The total torque M of the turbine is computed including summing, from rmin to rmax, (C xi  sin β i +C yi  cos β i ).

CROSS-REFERENCE TO RELATED APPLICATIONS

This application claims the benefit under 35 U.S.C. §119(e) of thefiling date of U.S. Provisional Patent Application No. 61/800,542, filedon Mar. 15, 2013, entitled “Enhanced Disk Actuator Modeling,” theentirety of which is herein incorporated by reference.

BACKGROUND

This specification relates to computer-implemented modeling of gas andliquid flows in turbines.

Computer simulation of gas and liquid flows in axial-flow turbo machinesis a sophisticated class of applied computational fluid dynamicsproblems due to combination of geometrical complexity of thecomputational domain and the large impact of additional physicaleffects, such as flow turbulence and separation, blade tip vertices, andwake instability. High fidelity discrete models based on 3D unsteadyNavier-Stokes equations result in computationally demanding andexpensive simulations.

Simpler models for such turbine systems can be constructed by a diskactuator model using 2D Navier-Stokes equations in which all bladewheels are modeled in an unducted turbine of zero thickness in anincompressible fluid of constant density.

SUMMARY

An Energy Extractor Disk Actuator model for simulating axial-flowpropellers, fans, turbine and turbo compressor impellers and bladewheels is suggested along with qualitative evaluation of its fidelity.The model is intended for preliminary cost-effective engineeringanalysis of turbo machines based on 2D (axially symmetric) Euler orNavier-Stokes equations. It uses the same approximate representation ofa blade wheel as an “equivalent” axially symmetric disk actuator asclassical Betz' model [1] and Glauert's Blade Element Momentum Theory[2]. However, the suggested model is capable of supporting asignificantly wider range of practical applications, since a) it doesnot use initial assumptions regarding magnitude of radial flow speed, b)it takes into account finite thickness of a blade wheel, and c) it isequally applicable to compressible flows, in contrast to the Betz' andGlauert's approaches.

In general, one innovative aspect of the subject matter described inthis specification can be embodied in methods that include the actionsof computing a total torque of a turbine from lift and drag coefficientsof the turbine, including obtaining, along multiple points of a blade ofa turbine from a minimum radius rmin of the blade to a maximum radiusrmax of the blade, lift coefficients C_(yi) and drag coefficientsC_(xi); obtaining, at the multiple points of the blade from rmin tormax, corresponding components of an upstream fluid flow velocity vectoru_(h,Ri) and u_(φ,Ri) and components of a downstream fluid flow velocityu_(h,Li) and u_(φ,Li); computing averaged directions β_(i), of theupstream and downstream fluid flow velocity vectors using the componentsof the upstream fluid flow velocity vector u_(h,Ri) and u_(φ,Ri) and thecomponents of the downstream fluid flow velocity u_(h,Li) and u_(φ,Li);and computing the total torque M of the turbine including summing, fromrmin to rmax, (C_(xi) sin β_(i)+C_(yi) cos β_(i)). Other embodiments ofthis aspect include corresponding computer systems, apparatus, andcomputer programs recorded on one or more computer storage devices, eachconfigured to perform the actions of the methods. A system of one ormore computers can be configured to perform particular operations oractions by virtue of having software, firmware, hardware, or acombination of them installed on the system that in operation causes orcause the system to perform the actions. One or more computer programscan be configured to perform particular operations or actions by virtueof including instructions that, when executed by data processingapparatus, cause the apparatus to perform the actions.

The foregoing and other embodiments can each optionally include one ormore of the following features, alone or in combination. The actionsinclude obtaining, from rmin to rmax, averaged dynamic pressures q_(i);using corresponding components of the upstream fluid flow velocityvector u_(h,Ri) and u_(φ,Ri) and corresponding components of thedownstream fluid flow velocity u_(h,Li) and u_(φ,Li), wherein the totaltorque M of the turbine is further based on the averaged dynamicpressures q_(i). computing the total torque M of the turbine comprisessumming, from rmin to rmax, q_(i) (C_(xi) sin β_(i)+C_(yi) cos β_(i)).The actions include obtaining a number of blades N of the turbine;obtaining, from rmin to rmax, chord lengths L_(i) of the blade; andcomputing, from rmin to rmax, solidity factors σ_(i) according toL_(i)×N/2π·r, wherein r is the length of the radius from rmin to rmax,wherein the total torque M of the turbine is further based on thesolidity factors σ_(i). Computing the total torque M of the turbinecomprises summing, from rmin to rmax, σ_(i)·(C_(xi) sin β_(i)+C_(yi) cosβ_(i)). The actions include computing the output mechanical power of theturbine P by multiplying the total torque M by turbine rotational speedΩ. The actions include computing, from rmin to rmax, enthalpy jumpvalues from an upstream enthalpy value H_(Ri) to a downstream enthalpyvalue H_(Li); and computing the total average power of the turbineP_(avg) including summing, from rmin to rmax, the enthalpy jump values.P_(avg) is given by:

P_(avg) = 2Π∫_(rm i n)^(rmax)r ⋅ ρ ⋅ u_(hi)(H_(R) − H_(L)) ⋅ 𝕕r,where r is the radius of the turbine blade. The actions includecomputing the coefficient of hydrodynamic efficiency η according to:

${\eta = \frac{M \times \Omega}{P_{avr}}},$wherein Ω is the turbine rotational speed. The enthalpy jump values foran incompressible fluid are given by:

${{H_{Ri} - H_{Li}} = {\frac{p_{Ri} - p_{Li}}{\rho} + \frac{u_{\varphi,{Ri}}^{2} - u_{\varphi,{Li}}^{2}}{2}}},$wherein, from rmin to rmax, p_(Ri) are upstream pressure values, p_(Li)are downstream pressure values, and p is the density of the fluid. Theenthalpy jump values for a compressible fluid are given by:

${{H_{Ri} - H_{Li}} = {\left\lbrack {E_{{int},R} + \frac{p_{Ri}}{\rho_{Ri}} + \frac{u_{h,{Ri}}^{2} + u_{\varphi,{Ri}}^{2}}{2}} \right\rbrack - \left\lbrack {E_{{int},L} + \frac{p_{Li}}{\rho_{Li}} + \frac{u_{h,{Li}}^{2} + u_{\varphi,{Li}}^{2}}{2}} \right\rbrack}},$wherein E_(int) represents the internal energy of the fluid per unitmass and depends on pressure p and temperature of the fluid T, andwherein P_(Ri) are upstream pressure values, p_(Li) are downstreampressure values, ρ_(Ri) are upstream density values, ρ_(Li) aredownstream density values of the fluid.

In general, another innovative aspect of the subject matter described inthis specification can be embodied in methods that include the actionsof computing a total torque of a turbine from lift and drag coefficientsof the turbine, including obtaining, along multiple points of a blade ofa turbine from a minimum radius rmin of the blade to a maximum radiusrmax of the blade, lift coefficients C_(yi) and drag coefficientsC_(xi); obtaining, at the multiple points of the blade from rmin tormax, corresponding components of an upstream fluid flow velocity vectoru_(hi) and u_(φI) and components of a downstream fluid flow velocityu_(hi) and u_(φI); computing averaged directions β_(i) of the upstreamand downstream fluid flow velocity vectors using the components of theupstream fluid flow velocity vector u_(hi) and u_(φI) and the componentsof the downstream fluid flow velocity u_(hi) and u_(φI); and computingthe total drag X of the turbine including summing, from rmin to rmax,(C_(yi) sin β_(i)−C_(xi) cos β_(i)). Other embodiments of this aspectinclude corresponding computer systems, apparatus, and computer programsrecorded on one or more computer storage devices, each configured toperform the actions of the methods.

The foregoing and other embodiments can each optionally include one ormore of the following features, alone or in combination. The actionsinclude obtaining, from rmin to rmax, averaged dynamic pressures q_(i);using corresponding components of the upstream fluid flow velocityvector u_(hi) and u_(φI) and corresponding components of the downstreamfluid flow velocity u_(hi) and u_(φI), wherein the total torque M of theturbine is further based on the averaged dynamic pressures q_(i).Computing the total drag X of the turbine comprises summing, from rminto rmax, q_(i)·(C_(yi) sin β_(i)−C_(xi) cos β_(i)). The actions includeobtaining a number of blades N of the turbine; obtaining, from rmin tormax, chord lengths L_(i) of the blade; and computing, from rmin tormax, solidity factors σ_(i) according to L_(i)×N/2π·r, wherein r is thelength of the radius from rmin to rmax, wherein the total torque M ofthe turbine is further based on the solidity factors σ_(i). Computingthe total drag X of the turbine comprises summing, from rmin to rmax,σ_(i)·(C_(yi) sin β_(i)−C_(xi) cos β_(i)).

Particular embodiments of the subject matter described in thisspecification can be implemented so as to realize one or more of thefollowing advantages. A system can compute more accurate parameters of aturbine without the computational complexity of a full 3D model. Thesystem can use 2D equations that are more accurate than a classical diskactuator model. The system can compute turbine parameters for ducted orunducted turbines, multistage turbines, turbines with stators orrotating blades, compressible or incompressible fluids, and turbineshaving a finite thickness.

The details of one or more embodiments of the subject matter of thisspecification are set forth in the accompanying drawings and thedescription below. Other features, aspects, and advantages of thesubject matter will become apparent from the description, the drawings,and the claims.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 illustrates the lateral section of a wide angle ducted turbinewith slot for boundary layer control.

FIG. 2 illustrates the integration area and contour in a lateral sectionof a disk actuator.

FIG. 3 is a flow chart of an example process for computing torque, drag,and power using jump conditions.

FIG. 4 illustrates the cross-sectional blade element at radius r andφ=π/2.

FIG. 5 illustrates a transformed integration area and contour in alateral section of disk actuator.

DETAILED DESCRIPTION

Computer simulation of gas and liquid flows in axial-flow turbo machinesconstitutes one of the most sophisticated classes of appliedcomputational fluid dynamics (CFD) problems due to combination ofgeometrical complexity of computational domain, essentially unsteadyflow pattern, and large impact of additional physical effects, such asflow turbulence and separation, blade tip vertices, wake instability,etc. So, high fidelity discrete models based on 3D unsteadyNavier-Stokes equations result in extremely computationally demandingand expensive simulation (see, for example, [3-5]). At preliminarystages of turbine design, for geometrically optimizing blade shapes andother components of turbine assembly it is necessary to run multipleproblem variants, so such precise models are unusable because ofenormous amount of required computer resources. In this case, asimplified (approximate) flow model could offer an appropriatecost-effective modeling tool for preliminary design.

Approximate models for simulating turbines are usually constructed onthe basis of various versions of the Energy Extractor Disk Actuatorprinciple originally suggested by A. Betz [1], and later extended by H.Glauert [2]. The Blade Element Momentum Theory (BEM) intended forapproximately representing hydrodynamic impact of an “equivalent” diskactuator on flow in a turbine in terms of hydrodynamic properties ofblade sections, had been also developed by the latter researcher.Present-day numerical techniques aimed at simplified design of turbinesuse the classic disk actuator principle and BEM in combination withanalytic or semi-analytic approximations of upstream and downstreamflows (see, for example, [6]), or numerically simulating upstream anddownstream flows using 2D Navier-Stokes equations (see, for example,[7]). In a number of most recent researches the disk actuator principleor its further extension for fully 3D flows, the Actuator Line Model(ALM), is used for modeling blade wheels in combination with full 3DReynolds-averaged Navier-Stokes equations (see, for example, [5,7]).However, as long as such hybrid approaches require solving a full 3D CFDproblem with turbulence and/or large eddy models, i.e. restore originalcomputational complexity, they cannot be considered as cost-effectivesolutions.

A more general axially symmetric disk actuator model is suggested below.Due to taking into account extra physical effects, such as finite radialflow speed, finite thickness of disk actuator, and fluidcompressibility, it provides significantly wider range of practicalapplications as compared with previously developed similar models.

The system of conservation laws describing fluid motion in a generalizedcurvilinear time-dependent coordinates (ξ, η, ζ), i.e. Navier-Stokesequations, can be written in the following divergent form:

$\begin{matrix}{{{\frac{\partial Q}{\partial t} + \frac{\partial\left( {E - E_{v}} \right)}{\partial\xi} + \frac{\partial\left( {G - G_{v}} \right)}{\partial\eta} + \frac{\partial\left( {F - F_{v}} \right)}{\partial ϛ}} = {S - S_{v}}},} & (1)\end{matrix}$

where Q is vector of mass, momentum, and energy densities; E, G, and Fare vectors of convective fluxes in respective directions ξ, η, and ζ;E_(v), G_(v), and F_(v) are vectors of dissipative fluxes; S and S_(v)are vectors of densities of respective convective and dissipativesources. Explicit expressions for those vectors depend on both physicalproperties of the fluid and selected coordinate system (ξ, η, ζ). If alldissipative fluxes and sources are zero E_(v)=0, G_(v)=0, F_(v)=0, andS_(v)=0, then Navier-Stokes equations (1) reduce to the Euler ones.

When using a differential turbulence model, system (1) should beextended with auxiliary equations describing evolution of localparameters of turbulence, such as average kinetic energy of turbulentmotion k and its dissipation rate ε in the (k−ε) model. Also, whenmodeling a flow of non-equilibrium chemically reacting gas mixture,system (1) should be extended with extra equations describing transitionand diffusion of mass concentrations of chemical components with sourceterms defining rates of production of the components in chemicalreactions. However, in our considerations below it is supposed thatturbulence and non-equilibrium chemical reactions do not createsignificant impact on free stream flow structure outside of boundarylayers, so the mentioned additional equations are not required.

Let's now consider typical problem of approximately modeling anaxial-flow turbo machine. An example of such problem, ducted hydroturbine with boundary layer control, is illustrated in FIG. 1 below. Aslong as simplified engineering model implies steady-state axiallysymmetric flow field, the latter can be compactly represented incylindrical coordinate system (h, r, φ) whose cylindrical axis coincideswith the turbine axis (see FIG. 1). In this case ξ=h, η=r, ζ=φ, ∂Q/∂t=0,∂(F−Fv)/∂ζ=0, and Navier-Stokes equations (1) take the following form:

$\begin{matrix}{{\frac{\partial\left( {E - E_{v}} \right)}{\partial h} + \frac{\partial\left( {G - G_{v}} \right)}{\partial r}} = {S - S_{v}}} & (2)\end{matrix}$

At high Reynolds numbers dissipative terms E_(v), G_(v), S_(v) arenegligibly small and can be omitted in free stream areas, i.e. outsideof boundary layers.

When constructing energy extractor disk actuator model, all blade wheelsof a turbo machine are replaced with azimuth-uniform disks of finite orzero thicknesses extracting the same average amounts of momentum andenergy as real blade wheels. For example, ducted hydro turbine in FIG. 1contains two such wheels: turbine impeller and guide vane wheel, so thattwo respective disk actuators should be used. For multi-staged turbinesand turbo compressors each rotor and stator blade wheel should bemodeled by a disk actuator, etc.

Therefore, each disk actuator approximately represents stream volumeswept by respective blade wheel and works as source or sink of momentumand energy. When using disk actuators in combination with axiallysymmetric Navier-Stokes equations (2), it is supposed that equations (2)are responsible for description of the flow field everywhere outside ofswept volumes, while disk actuators provide “jump conditions” connectingflow field parameters at upstream (“L”) and downstream (“R”) sides ofrespective blade wheels. Derivation of the jump conditions is presentedbelow.

Vectors of convective fluxes E, G and sources S in axially symmetricNavier-Stokes equations (2) can be expressed in terms of flow parametersas follows:

$\begin{matrix}{{E = {\begin{matrix}\begin{matrix}\begin{matrix}\begin{matrix}{r\;\rho\; u_{h}} \\{r\left( {{\rho\; u_{h}^{2}} + p} \right)}\end{matrix} \\{r\;\rho\; u_{h}u_{r}}\end{matrix} \\{r^{2}\rho\; u_{h}u_{\varphi}}\end{matrix} \\{r\;\rho\; u_{h}H^{{(*})}}\end{matrix}}},{G = {\begin{matrix}\begin{matrix}\begin{matrix}\begin{matrix}{r\;\rho\; u_{r}} \\{r\;\rho\; u_{h}u_{r}}\end{matrix} \\{r\left( {{\rho\; u_{r}^{2}} + p} \right)}\end{matrix} \\{r^{2}\rho\; u_{r}u_{\varphi}}\end{matrix} \\{r\;\rho\; u_{r}H^{{(*})}}\end{matrix}}},{S = {\begin{matrix}\begin{matrix}\begin{matrix}\begin{matrix}0 \\{rf}_{h}\end{matrix} \\{p + {\rho\; u_{\varphi}^{2}} + {rf}_{r}}\end{matrix} \\{r^{2}f_{\varphi}}\end{matrix} \\{re}\end{matrix}}}} & (3)\end{matrix}$

where uh, ur, and uφ are components of velocity vector in respectivedirections h, r, and φ, ρ is fluid density, p is static pressure, H(*)is total specific enthalpy of fluid per unit mass. Source terms fh, fr,and fφ represent volume densities of components of external forces, ande is volume density of energy sources.

Although expressions (3) are valid for both compressible andincompressible fluids, standard definitions of the total specificenthalpy H(*) and physical meaning of the energy equation are different.For incompressible fluids enthalpy H* is usually introduced totalmechanical energy, i.e. sum of potential and kinetic energies, per unitmassH*=p/ρ+u ²/2=p/ρ+(u _(h) ² +u _(r) ² +u _(φ) ²)/2.  (4)

On the other hand, definition of enthalpy H for compressible fluidsincludes additional term representing internal (thermal) energy per unitmassH=E _(int) +p/ρ+u ²/2=E _(int) +p/ρ+(u _(h) ² +u _(r) ² +u _(φ)²)/2,  (5)

where E_(int)(p, T) is internal energy of the fluid per unit mass, and Tis fluid temperature. For example, in case of a thermally perfect gasE_(int)=C_(v)T/μ and (5) reduces toH=C _(v) T/μ+p/ρ+u ²/2=C _(p) T/μ+(u _(h) ² +u _(r) ² +u _(φ) ²)/2,  (6)

where C_(v) and C_(p) are specific heat capacities at constant volumeand pressure, respectively, and μ is molar mass of the gas. System ofNavier-Stokes equations (1) for a compressible gas includes densityρ(t,ξ,η,ζ) as an additional unknown function as compared with system (1)for incompressible fluid, so equation of gas state ρ=ρ(p, T) has to beused for closing the problem.

As long as expression (4) does not include internal energy Eint, it doesnot represent a conservative physical value. Therefore, if energyequations for an incompressible fluid is written in terms of thefunction (4), it does not reflect an independent conservation law,because H* may dissipate due to the impact of viscous friction and heattransfer (Sv≠0) at finite Reynolds numbers. In fact, it can be shownthat the considered equation directly follows from the mass and momentumconservation laws. On the other hand, energy equation for a compressiblegas written in terms of the function (5) or (6) does represent anindependent conservation law defining temperature distribution.

Vector equation (2) can be re-written in the scalar form(∂E _(k) /∂h)+(∂G _(k) /∂r)=S _(k),  (7)where Ek, Gk and Sk denote k-th scalar components of the vectors(E−E_(v)), (G−G_(v)), and (S−S_(v)), respectively. Index k in (7) variesfrom 1 to 5, thus representing continuity equation (k=1), h-momentumequation (k=2), r-momentum equation (k=3), φ-momentum equation (k=4),and energy equation (k=5). Let's now select an infinitely narrowintegration area A in plane (h, r), so that its lower and upper boundscross lateral section of disk actuator and are located at radii r andr+Δr, respectively, while its left (h=h_(L)(r)) and right (h=h_(R)(r))bounds are located in the areas of free stream closely adjoining leftand right sides of the disk, see FIG. 2. In our considerations below,flow parameters at the left and right bounds of the area A will bemarked with the indices “L” and “R”, respectively.

Integrating scalar equations (7) over area A and applying the 1-stGreen's integral formula, yields

$\begin{matrix}{{{\oint\limits_{C}\left( {{E_{k} \cdot {\mathbb{d}r}} - {G_{k} \cdot {\mathbb{d}h}}} \right)} = {\underset{A}{\int\int}{S_{k} \cdot {\mathbb{d}r} \cdot {\mathbb{d}h}}}},} & (8)\end{matrix}$

where C is the contour bounding area of integration A. Since height Δrof the area A is supposed to be infinitely small (Δr→0), the contour andarea integrals in equation (8) can be approximated with an arbitraryprecision as follows:

$\mspace{20mu}{{{\oint\limits_{C}{E_{k} \cdot {\mathbb{d}r}}} = {\Delta\;{r \cdot \left\{ {\left\lbrack {E_{k}(r)} \right\rbrack_{R} - \left\lbrack {E_{k}(r)} \right\rbrack_{L}} \right\}}}},{{\oint\limits_{C}{G_{k} \cdot {\mathbb{d}h}}} = {{\int_{h_{L}{(r)}}^{h_{R}{(r)}}{\left\lbrack {{G_{k}\left( {h,r} \right)} - {G_{k}\left( {h,{r + {\Delta\; r}}} \right)}} \right\rbrack \cdot {\mathbb{d}h}}} = {{{- \Delta}\;{r \cdot {\int_{h_{L}{(r)}}^{h_{R}{(r)}}{\left\lbrack \frac{\partial{G_{k}(r)}}{\partial r} \right\rbrack \cdot {\mathbb{d}h}}}}} = {{- \Delta}\;{{rD}\left\lbrack \frac{\partial{G_{k}(r)}}{\partial r} \right\rbrack}_{avr}}}}},\mspace{20mu}{{\underset{A}{\int\int}{S_{k} \cdot {\mathbb{d}r} \cdot {\mathbb{d}r}}} = {{\Delta\;{r \cdot {S_{k}\left( {h,r} \right)} \cdot {\mathbb{d}h}}} = {\Delta\;{{rD}\left\lbrack {S_{k}(r)} \right\rbrack}_{avr}}}},}$

where D=D(r)=h_(R)(r)−h_(L)(r) is local thickness of disk actuator,index “avr” means respective average values on the segmenth_(L)(r)≦h≦h_(R)(r), and indices “L” and “R” denote flow parameter atthe “left” and “right” points (h, r) closely adjoining upstream anddownstream sides of actuator, respectively, as indicated in FIG. 2.After substituting approximations above in equation (8), the latterreduces to[E _(k)(r)]_(R) −[E _(k)(r)]_(L) =D[S _(k)(r)−∂G _(k)(r)/∂r]_(avr).  (9)

Equation (9) can now be directly used for deriving jump conditions forflow parameters.

Since left and right bounds of A are located in free stream areas, i.e.outside of boundary layers, contribution of dissipative terms to theaxial fluxes [E_(k)(r)]_(R) and [E_(k)(r)]_(L) is negligibly small athigh Reynolds numbers, see note 1 above. On the other hand, impact ofthe dissipative effects on the average values [S_(k)(r)]_(avr) and[∂G_(k)(r)/∂r]_(avr) inside swept areas of blade wheels can besignificant, but when applying disk actuator model, it is assumed thatall such effects are represented by properly defined “effective” dragand lift coefficients C_(D) and C_(L) of blade elements defining sourceterms f_(h), f_(r), f_(φ), and e in (3), see equations 28-49 below.Therefore, explicitly including dissipative terms in expressions forE_(k), G_(k) and S_(k) is not necessary.

For evaluating accuracy of approximate formulas the standard “big O”notation will be used below. Relationship a=O(b) between two physicalvalues a and b should be interpreted as “a has the same order ofmagnitude as b, or less”, so that ratio |a/b| is always finite. Also,two dimensionless parameters will be used:ε=u _(r) /u _(h) , δ=D/r.  (10)

which may be small or finite depending on a particular flow pattern andconsidered point of flow field. As long as all previous considerationswere made in terms of dimensional variables, dimensional velocity scaleu∞ and density scale ρ_(∞) are also required for comparing orders ofmagnitudes of physical values. For example, when analyzing an open windor hydro turbine, velocity and density of the far upstream flow can beaccepted as such scales. Generally, it is supposed that scales u∞, ρ∞are selected so thatu _(h) =u _(∞) ·O(1), ρ=ρ_(∞) ·O(1)  (11)

everywhere inside swept volume of a blade wheel. Relationship (12)results in the following initial evaluations:u _(r) =u _(∞) ·O(ε), u _(φ) =u _(∞) ·O(1),Δp=ρ _(∞) u _(∞) ² ·O(1), ΔH=u _(∞) ² ·O(1).  (12)

In flows of an incompressible fluid static pressure p is defined withaccuracy of an arbitrary spatially constant additive function p₀(t), sothat grad [p(t,x,y,z)+p₀(t)]=grad p(t,x,y,z). In such cases only spatialvariations of pressure Δp, or grad(p) make physical sense, but not itsabsolute values.

For k=1 scalar components of vectors (3) areE₁=rρu_(h), G₁=rρu_(r), S₁=0.

So, the 1st of averaged equations (9) reflecting mass balance acrossdisk actuator isr[(ρu _(h))_(R)−(ρu _(h))_(L) ]=−D[∂(rρu _(r))/∂r] _(avr).

As long as ∂(rρu_(r))/∂r=ρ_(∞)u_(∞)·O(ε), this finally results inestimation(ρu _(h))_(R)−(ρu _(h))_(L)=ρ_(∞) u _(∞) ·O(ε·δ).  (13)

In case of incompressible fluid ρ=ρ_(∞)=const and estimation (13)reduces to(u _(h))_(R)−(u _(h))_(L) =u _(∞) ·O(ε·δ).  (14)

For k=2 scalar components of vectors (3) areE ₂ =r(ρu _(h) ² +p), G ₂ =rρu _(h) u _(r) , S ₂ =rf _(h),

where dissipative terms (E_(v))₂, (G_(v))₂, and (S_(v))₂ are omitted inaccordance with note 3 above. So, the 2nd of averaged equations (9)reflecting h-momentum balance across disk actuator isr[(ρu _(h) ² +p)_(R)−(ρu _(h) ² +p)_(L) ]=D[rf _(h)−∂(rρu _(h) u_(r))/∂r] _(avr).As long as ∂(rρu_(h)u_(r))/∂r=ρ_(∞)u_(∞) ²·O(ε), this finally results inestimation(ρu _(h) ² +p)_(R)−(ρu _(h) ² +p)_(L) =D·(f _(h))_(avr)+ρ_(∞) u _(∞) ²·O(ε·δ),

and taking into account (13),(u _(h))_(R)−(u _(h))_(L)+(p _(R) −p _(L))/(ρu _(h))=D·(f_(h))_(avr)/(ρu _(h))+u _(∞) ·O(ε·δ),  (15)

In case of incompressible fluid combination of estimations (14) and (15)yieldsp _(R) −p _(L) =D·(f _(h))_(avr)+ρ_(∞) u _(∞) ² ·O(ε·δ).  (16)

Note that average volume density of axial forces (f_(h))_(avr) in (15)and (16) may not be finite, since for a finite total impact of diskactuator on h-momentum balance D·(f_(h))_(avr)=O(1), hence|(f_(h))_(avr)|→∞ at D→0.

For k=3 scalar components of vectors (3) areE ₃ =rρu _(h) u _(r) , G ₃ =r(ρu _(r) ² +p), S ₃ =p+ρu _(φ) ² +rf _(r),

where dissipative terms (E_(v))₃, (G_(v))₃, and (S_(v))₃ are omitted inaccordance with note 3 above. So, the 3rd of averaged equations (9)reflecting r-momentum balance across disk actuator is

$\begin{matrix}{{r\left\lbrack {\left( {\rho\; u_{h}u_{r}} \right)_{R} - \left( {\rho\; u_{h}u_{r}} \right)_{L}} \right\rbrack} = {{D\left\{ {p + {\rho\; u_{\varphi}^{2}} + {rf}_{r} - \frac{\partial\left\lbrack {r\left( {{\rho\; u_{r}^{2}} + p} \right)} \right\rbrack}{\partial r}} \right\}_{avr}} =}} \\{{= {D\left\lbrack {{\rho\; u_{\varphi}^{2}} + {rf}_{r} - {r \cdot \frac{\partial p}{\partial r}} - \frac{\partial\left( {r\;\rho\; u_{r}^{2}} \right)}{\partial r}} \right\rbrack}_{avr}},}\end{matrix}$

and after dividing by r we get(ρu _(h) u _(r))_(R)−(ρu _(h) u _(r))_(L) =δ[ρu _(φ) ² +rf _(r)−r·∂p/∂r−∂(rρu _(r) ²)/∂r] _(avr).  (17)

As long average volume density of radial forces (f_(r))_(avr) is lessthen (f_(h))_(avr) and (fφ)_(avr) by its order of magnitude, and aninfinitely thin blade wheel should provide zero impact on r-momentum ofthe flow, i.e. D·(f_(h))_(avr)→0 when D→0, it would be reasonable toassume that (f_(r))_(avr) remains finite independently of D, i.e.(f_(r))_(avr)=O(1). If this assumption is valid, then the right handside term in square brackets is finite as well:[ρu _(φ) ² +rf _(r) −r·∂p/∂r−∂(rρu _(r) ²)/∂r] _(avr) =O(1)

everywhere in the area of integration A, and equation (17) immediatelyresults in estimation(ρu _(h) u _(r))_(R)−(ρu _(h) u _(r))_(L)=ρ_(∞) u _(∞) ² ·O(δ).  (18)On the other hand, u_(r)=u_(∞)·O(ε) by definition (10) of parameter ε,so that(ρu _(h) u _(r))_(R)−(ρu _(h) u _(r))_(L)=ρ_(∞) u _(∞) ² ·O(ε·δ).  (19)

Combination of (18) and (19) yields(ρu _(h) u _(r))_(R)−(ρu _(h) u _(r))_(L)=ρ_(∞) u _(∞) ² ·O(ε·δ),

and taking into account jump condition (13), we finally get(u _(r))_(R)−(u _(r))_(L) =u _(∞) ·O(ε·δ),  (20)

for both compressible and incompressible fluids.

For k=4 scalar components of vectors (3) areE₄=r²ρu_(h)u_(φ), G₄=r²ρu_(r)u_(φ), S₄=r²f_(φ),

where dissipative terms (E_(v))₄, (G_(v))₄, and (S_(v))₄ are omitted inaccordance with note 3 above. So, the 4th of averaged equations (9)reflecting φ-momentum balance across disk actuator isr ²[(ρu _(h) u _(φ))_(R)−(ρu _(h) u _(φ))_(L) ]=D[r ² f _(φ)−∂(r ² ρu_(r) u _(φ))/∂r] _(avr).

As long as [∂(r²ρu_(r)u_(φ))/∂r]_(avr)=rρ_(∞)u_(∞) ²·O(ε), this resultsin estimation(ρu _(h) u _(φ))_(R)−(ρu _(h) u _(φ))_(L) =D·(f _(φ))_(avr)+ρ_(∞) u _(∞)² ·O(ε·δ),

and taking into account jump condition (13), we finally get(u _(φ))_(R)−(u _(φ))_(L) =D·(f _(φ))_(avr)/(ρu _(h))+u _(∞)·O(ε·δ)  (21)

for both compressible and incompressible fluids. Note that averagevolume density of azimuth forces (f_(φ))_(avr) in (21) may not befinite, since for a finite total impact of disk actuator on φ-momentumbalance D·(f_(φ))_(avr)=O(1), hence |(f_(φ))_(avr)|→∞ at D→0.

For k=5 scalar components of vectors (3) areE₅=rρu_(h)H^((*)), G₅=rρu_(r)H^((*)), S₅=re,

where dissipative terms (E_(v))₅, (G_(v))₅, and (S_(v))₅ are omitted inaccordance with note 3 above. So, the 5th of averaged equations (9)reflecting energy balance across disk actuator isr[(ρu _(h) H ^((*)))_(R)−(ρu _(h) H ^((*)))_(L) ]=D[re−∂(rρu _(r) H^((*)))/∂r] _(avr).As long as ∂(rρu_(r)H^((*)))/∂r=ρ_(∞)u_(∞) ³·O(ε), this results inestimation(ρu _(h) H ^((*)))_(R)−(ρu _(h) H ^((*)))_(L) =D·e _(avr)+ρ_(∞) u _(∞) ³·O(ε·δ),

and taking into account jump condition (13), we finally getH ^((*)) _(R) −H ^((*)) _(L) =D·e _(avr)/(ρu _(h))+u _(∞) ²·O(ε·δ)  (22)

for both compressible and incompressible fluids. Note that averagevolume density of energy sources e_(avr) in (22) may not be finite,since for a finite total impact of disk actuator on energy balanceD·e_(avr)=O(1), hence |e_(avr)|→∞ at D→0.

FIG. 3 is a flow chart of an example process for computing torque, drag,and power using jump conditions. The process will be described as beingperformed by an appropriately programmed system of one or morecomputers.

The system obtains mass, momentum, and energy jump conditions (310).After summarizing estimations (13), (15), (20), (21), and (22), we canconclude that the approximate jump conditions connecting flow parameterson the upstream and downstream sides of disk actuator

$\begin{matrix}\left\{ \begin{matrix}\begin{matrix}\begin{matrix}\begin{matrix}{{\left( {\rho\; u_{h}} \right)_{R} = \left( {\rho\; u_{h}} \right)_{L}},} \\{{{\left( u_{H} \right)_{R} - \left( u_{h} \right)_{L} + {\left( {p_{R} - p_{L}} \right)/\left( {\rho\; u_{h}} \right)}} = {D \cdot {\left( f_{h} \right)_{avr}/\left( {\rho\; u_{h}} \right)}}},}\end{matrix} \\{{\left( u_{r} \right)_{R} = \left( u_{r} \right)_{L}},}\end{matrix} \\{{{\left( u_{\varphi} \right)_{R} - \left( u_{\varphi} \right)_{L}} = {D \cdot {\left( f_{\varphi} \right)_{avr}/\left( {\rho\; u_{h}} \right)}}},}\end{matrix} \\{{H_{R}^{{(*})} - H_{L}^{{(*})}} = {D \cdot {e_{avr}/\left( {\rho\; u_{h}} \right)}}}\end{matrix} \right. & \begin{matrix}\begin{matrix}\begin{matrix}\begin{matrix}(23) \\(24)\end{matrix} \\(25)\end{matrix} \\(26)\end{matrix} \\(27)\end{matrix}\end{matrix}$

are all valid with relative accuracy O(ε·δ). Therefore, fidelity of themodel is determined by the dimensionless scale factor (ε·δ), rather thanscale of radial velocity ε, i.e. initial assumption furl uh, is notrequired. Instead, the assumption (ε·δ)<<1 should be used as a fidelitycriterion, which is satisfied for virtually all practical cases.

In fact, in all axial-flow turbo machines root parts of blades (largeδ=D/r) do not create significant impact on momentum and energy of theflow because of their small swept areas, tangential speed, and armlength of applied hydrodynamic loads. On the other hand, middle and tipparts of blades (small δ=D/r), producing major contribution to momentumand energy re-distribution, satisfy to the condition (ε·δ)<<1,especially as their pitch angles are typically close to 90°.

The average h- and φ-momentum sources (f_(h))_(avr), (f_(φ))_(avr) inright hand sides of jump conditions (24), (26) can be evaluated it termshydrodynamic properties of 2D blade sections using classical Glauert'sBlade Element Momentum Model (BEM), described, for example, indissertation [7].

For convenience of further considerations, an auxiliary Cartesiancoordinate system (x, y, z) will be used below in addition to thepreviously introduced cylindrical coordinates (h, r, φ):

$\begin{matrix}\left\{ \begin{matrix}{{x = h},} \\{{y = {{r \cdot \cos}\;\varphi}},} \\{z = {{r \cdot \sin}\;{\varphi.}}}\end{matrix} \right. & (28)\end{matrix}$

Cartesian components of flow velocity u, v, w, in respective directionsx, y, z can be expressed in terms of u_(n), u_(r), u_(φ) as

$\begin{matrix}\left\{ \begin{matrix}{u = u_{h}} \\{v = {{{u_{r} \cdot \cos}\;\varphi} - {{u_{\varphi} \cdot \sin}\;\varphi}}} \\{w = {{{u_{r} \cdot \sin}\;\varphi} + {{u_{\varphi} \cdot \cos}\;\varphi}}}\end{matrix} \right. & (29)\end{matrix}$

Let's now consider an infinitely thin blade element formed byintersection of blade with a cylinder of given radius r. Such elementlocated in azimuth position φ=π/2 is shown in FIG. 4, where all upstreamand downstream flow parameters, as well as blade section geometry,depend on radius r. Note that u=u_(h), v=−u_(φ), and w=u_(r) at φ=π/2.

In the coordinate system rotating with a constant speed Ω, where Ω isrotational speed of the blade wheel, Cartesian velocity components are(u, v−Ωr, w). Let's introduce additional local Cartesian coordinates(x′, y′) rotated by angle β_(avr) with respect to initial coordinates(x, y) in tangential plane of the blade element:

$\begin{matrix}{{x^{\prime} = {{x\;\cos\;\beta_{avr}} + {y\;\sin\;\beta_{avr}}}},{y^{\prime} = {{{- x}\;\sin\;\beta_{avr}} + {y\;\cos\;\beta_{avr}}}},{x = {{x^{\prime}\;\cos\;\beta_{avr}} - {y^{\prime}\;\sin\;\beta_{avr}}}},{y = {{x^{\prime}\;\sin\;\beta_{avr}} + {y^{\prime}\;\cos\;\beta_{avr}}}},{where}} & (30) \\\begin{matrix}{\beta_{avr} = {{{atan}\left\lbrack {\left( {v - {\Omega\; r}} \right)_{avr}/u_{avr}} \right\rbrack} =}} \\{= {{- {atan}}\left\{ \frac{\left\lbrack {\left( u_{\varphi} \right)_{L} + \left( u_{\varphi} \right)_{R} + {2\Omega\; r}} \right\rbrack}{\left\lbrack {\left( u_{h} \right)_{L} + \left( u_{h} \right)_{R}} \right\rbrack} \right\}}}\end{matrix} & (31)\end{matrix}$

characterizes an averaged direction of flow velocity in rotatingcoordinate system. Formula (31) is specifically constructed to providesuch direction of total hydrodynamic load applied to the blade elementthat for an ideal hydrofoil having zero drag, amount of power withdrawnby the element is exactly zero in rotating coordinate system, i.e. thereis no dissipative energy losses. If the blade element produces positiveoutput power then Ωr>|v_(L)|, so that β_(avr)<0.

Components of hydrodynamic load applied to the blade element in thecoordinate system (x′, y′) areΔX′=[(ρU ²)_(avr)/2]C _(D) L·Δr=qC _(D) L·Δr,ΔY′=[(ρU ²)_(avr)/2]C _(L) L·Δr=qC _(L) L·Δr,  (32)

where Δr is radial thickness of the element (Δr→0), L=L(r) is its chord,C_(D)=C_(D)(α, r), C_(L)=C_(L)(α, r) are its drag and lift coefficients,respectively, and α=θ(r)+β_(avr)(r) is its angle of attack in the localcoordinate system (x′, y′). The average dynamic pressure q=(ρU²)_(avr)/2in (32) is defined as

$\begin{matrix}\begin{matrix}{q = {\left( {\rho\; U^{2}} \right)_{avr}/2}} \\{= {\frac{\left\{ {\rho\left\lbrack {u^{2} + \left( {v - {\Omega\; r}} \right)^{2}} \right\rbrack} \right\}_{avr}}{2} =}} \\{= {\frac{\begin{Bmatrix}{{\rho_{L}\left\lbrack {\left( u_{h} \right)_{L}^{2} + \left( {u_{\varphi} + {\Omega\; r}} \right)_{L}^{2}} \right\rbrack} +} \\{\rho_{R}\left\lbrack {\left( u_{h} \right)_{R}^{2} + \left( {u_{\varphi} + {\Omega\; r}} \right)_{R}^{2} + \left( {u_{\varphi} + {\Omega\; r}} \right)_{R}^{2}} \right\rbrack}\end{Bmatrix}}{4}.}}\end{matrix} & (33)\end{matrix}$

Radial velocity u_(r) is not included in right hand side of expression(33) as its contribution to hydrodynamic loads is supposed to benegligible, see equation 17.

Components of hydrodynamic load applied to the blade element in theinitial coordinate system (x, y) areΔX=ΔX′ cos β_(avr) −ΔY′ sin β_(avr) =q(C _(D) cos β_(avr) −C _(L) sinβ_(avr))L·Δr,  (34)ΔY=ΔX′ sin β_(avr) +ΔY′ cos β_(avr) =q(C _(D) sin β_(avr) +C _(L) cosβ_(avr))L·Δr.  (35)

Drag, torque, and mechanical power produced by the blade element are:ΔX, r·ΔY, and Ωr·ΔY, respectively, so that their total values can becomputed via integration from r=r_(min) to r=r_(max), and multiplicationresulting integrals by number of blades.

For the considered blade element above, averaged values (f_(h))_(avr),(f_(φ))_(avr) from (24,26) can be evaluated from the forces (34,35) asfollows:(f _(h))_(avr) =−ΔX/ΔV=−ΔX/(2πrD·Δr/N),  (36)(f _(φ))_(avr) =ΔY/ΔV=ΔY/(2πrD·Δr/N),  (37)

where N is number of impeller blades, so that each element of eachseparate blade provides momentum input in the azimuth range 2π/N, sorespective elementary swept volume is ΔV=2πrD·Δr/N. As long as(f_(h))_(avr), (f_(φ))_(avr) in (24), (26) represent flow reaction tothe hydrodynamic loads applied to the blade element, they should beassigned opposite signs. That is why in formula (36) sign “−” is set forΔX (as x- and h-directions coincide, ∂x/∂h>0), while in formula (37)sign “+” is set for ΔY (as y- and φ-directions are opposite, ∂y/∂φ<0 atφ=π/2). Substitution of expressions (36), (37) in jump conditions (24),(26), respectively, and use of relationships (34), (35) for evaluatingratios ΔX/Δr, ΔY/Δr, results in representation of h- and (p-momentumjumps in terms of C_(D) and C_(L):

$\begin{matrix}\begin{matrix}{{{\rho\;{u_{h}\left\lbrack {\left( u_{h} \right)_{R} - \left( u_{h} \right)_{L}} \right\rbrack}} + \left( {p_{R} - p_{L}} \right)} = {{{- \left( {\Delta\;{X/\Delta}\; r} \right)}{N/\left( {2\pi\; r} \right)}} =}} \\{= {{- {q\begin{pmatrix}{{C_{D}\cos\;\beta_{avr}} -} \\{C_{L}\sin\;\beta_{avr}}\end{pmatrix}}}L\;{N/\left( {2\pi\; r} \right)}}} \\{{= {q\;{\sigma\left( {{C_{L}\sin\;\beta_{avr}} - {C_{D}\cos\;\beta_{avr}}} \right)}}},}\end{matrix} & (38) \\{\mspace{20mu}\begin{matrix}{{\rho\;{u_{h}\left\lbrack {\left( u_{\varphi} \right)_{R} - \left( u_{\varphi} \right)_{L}} \right\rbrack}} = {{\left( {\Delta\;{Y/\Delta}\; r} \right){N/\left( {2\pi\; r} \right)}} =}} \\{= {{q\left( {{C_{D}\sin\;\beta_{avr}} + {C_{L}\cos\;\beta_{avr}}} \right)}L\;{N/\left( {2\pi\; r} \right)}}} \\{{= {q\;{\sigma\left( {{C_{D}\sin\;\beta_{avr}} + {C_{L}\cos\;\beta_{avr}}} \right)}}},}\end{matrix}} & (39)\end{matrix}$

where σ=σ(r)=LN/(2πr) is solidity factor of blade wheel.

The system determines the upstream flow field parameters and thedownstream flow field parameters that satisfy the jump conditions (320),for example, by using a computational fluid dynamics software package.The system determines momentum jump values from the upstream flow fieldparameters and the downstream flow field parameters (330).

The system computes the drag, torque, and power produced by the bladewheel from the momentum jumps (340). Drag ΔX, torque ΔM and mechanicalpower ΔW produced by N blade elements areN·ΔX=−2πr{ρu _(h)[(u _(h))_(R)−(u _(h))_(L)]+(p _(R) −p _(L))}·Δr,N·ΔM=2πr ² ρu _(h)[(u _(φ))_(R)−(u _(φ))_(L) ]·Δr,N·ΔW=ΩN·ΔM,  (40)

and total drag X, torque M, and mechanical power W from the energy andproduced by the blade wheel are

$\begin{matrix}{{X = {{- 2}\pi{\int_{\underset{\_}{r_{m\; i\; n}}}^{\underset{\_}{r_{{ma}\; x}}}{r{\left\{ {{\rho\;{u_{h}\left\lbrack {\left( u_{h} \right)_{R} - \left( u_{h} \right)_{L}} \right\rbrack}} + \left( {p_{R} - p_{L}} \right)} \right\} \cdot {\mathbb{d}r}}}}}},{M = {2\pi{\int_{\underset{\_}{r_{m\; i\; n}}}^{\underset{\_}{r_{{ma}\; x}}}{r^{2}\rho\;{{u_{h}\left\lbrack {\left( u_{\varphi} \right)_{R} - \left( u_{\varphi} \right)_{L}} \right\rbrack} \cdot {\mathbb{d}r}}}}}},{W = {\Omega \cdot M}},} & (41)\end{matrix}$

where the integrand expressions are represented by (38) and (39).

Relationships (23), (25), (38), and (39) constitute a complete set ofjump conditions for equations of continuity and momentums. Jumpcondition (27) for energy equation should be considered separately forthe cases of incompressible and compressible fluid because ofdistinction of their physical meanings and different interpretations oftotal specific enthalpy H(*), see note 2 above.

If the energy equation, i.e. equation (7) for k=5, is written withouttaking into account internal (thermal) fluid energy, it describesdistribution of total specific enthalpy (4) representing only potential(p/ρ) and kinetic (u²/2) energies per unit mass. Such form of energyequation can be derived as a combination of all remaining equations, sothat density of energy sources e in (3) can be expressed in terms ofmomentum sources f_(h) and f_(φ). In this case jump of enthalpy (4) canbe expressed in terms of pressure and velocity components, without usingcondition (27), as follows:

$\begin{matrix}\begin{matrix}{{H_{R}^{*} - H_{L}^{*}} = {\left\lbrack {{p/\rho} + {\left( {u_{h}^{2} + u_{r}^{2} + u_{\varphi}^{2}} \right)/2}} \right\rbrack_{R} -}} \\{\left\lbrack {{p/\rho} + {\left( {u_{h}^{2} + u_{r}^{2} + u_{\varphi}^{2}} \right)/2}} \right\rbrack_{L} =} \\{= {{\left( {p_{R} - p_{L}} \right)/\rho} + {\left\lbrack {\left( u_{\varphi} \right)_{R}^{2} - \left( u_{\varphi} \right)_{L}^{2}} \right\rbrack/2}}}\end{matrix} & (42)\end{matrix}$

using previously computed p and u_(φ), and taking into account that(u_(r))_(R)=(u_(r))_(L) and (u_(h))_(R)=(u_(h))_(L) at constant ρ, see(23), (25). So, mechanical power contributed to the flow by N bladeelements located at radius r equalsN·ΔP*=2πrρu _(h)(H* _(R) −H* _(L))·Δr,  (43)

and total mechanical power P* contributed by the blade wheel is

$\begin{matrix}{P^{*} = {2\pi{\int_{r_{\min}}^{r_{\max}}{r\;\rho\;{{u_{h}\left( {H_{R}^{*} - H_{L}^{*}} \right)} \cdot \ {{\mathbb{d}r}.}}}}}} & (44)\end{matrix}$

Note that generally W≠−P* (see (41)), because the latter includes partof mechanical power dissipated into heat due to viscous friction.

Let's consider an ideal blade element having zero drag (C_(D)=0), andsuppose that input flow is irrotational ((u_(φ))_(L)=0). In this casejump conditions (38), (39) reduce to

$\begin{matrix}{\mspace{79mu}{{{{p_{R} - p_{L}} = {q\;\sigma\; C_{L}\;\sin\;\beta_{avr}}},\mspace{79mu}{\left( u_{\varphi} \right)_{R} = {q\;\sigma\; C_{L}\mspace{11mu}\cos\;{\beta_{avr}/\left( {\rho\; u_{h}} \right)}}},\mspace{20mu}{hence}}{{{p_{R} - p_{L}} = {{\rho\;{u_{h}\left( u_{\varphi} \right)}_{R}\;{\tan\left( \beta_{avr} \right)}} = {{{- {{\rho\left( u_{\varphi} \right)}_{R}\left\lbrack {\left( u_{\varphi} \right)_{R} + {2\Omega\; r}} \right\rbrack}}/2}=={{- \rho}\; r^{2}{\omega\left( {{\omega/2} + \Omega} \right)}}}}},}}} & (45)\end{matrix}$

where tan (β_(avr)) is evaluated using expression (31) at(u_(h))_(L)=(u_(h))_(R)=u_(h), and rotational speed ω_(φ)=(u_(φ))_(R)/ris introduced. Note that relationship (45) exactly coincides with theclassical Glauert's estimation of pressure jump across disk actuator. Inthe considered idealized case relationship (42) reduces toH* _(R) −H* _(L)=(p _(R) −p _(L))/ρ+r ²ω²/2=−r ²ωΩ,thus resulting in N·ΔW=Ω·2πr ²ρ_(h)[(u _(φ))_(R)−(u _(φ))_(L) ]·Δr=2πr ³ρu _(h) ωΩ·Δr,and N·ΔP*=−2πr ³ ρu _(h) ωΩ·Δr=−N·ΔW,

see (40) and (43). Therefore, W=−P* in the absence of dissipative losses(C_(D)=0).

When practically evaluating efficiency of open wind and water turbines,the impact of blade tip vortex structures on momentum and energybalances should also be taken into account. Contribution of the tipvortices in impeller wake to drag X and power P* of the turbine can beevaluated using one of standard tip correction techniques, which areminutely described in the publications specifically dedicated topractical implementations of disk actuator models, see for example [7].

If energy equation, i.e. equation (7) for k=5, is written with takinginto account internal (thermal) fluid energy, it describes distributionof total specific enthalpy (5) or (6) representing full energy of thefluid per unit mass. In this case density of energy sources e in (3)includes both power produced by momentum sources f_(h), f_(φ), and heatfluxes arising from non-zero gradients of temperature in a thermallyconductive media, so that jump condition (27) reflects full amount ofpower contributed to the flow. However, instead of directly evaluatingenergy source e and using condition (27), it is more convenient toderive an independent relationship of energy balance.

Jump of enthalpy (5) can be expressed in term of velocity components andother flow parameters, as follows:

$\begin{matrix}\begin{matrix}{{H_{R} - H_{L}} = {\left\lbrack {E_{int} + {p/\rho} + {\left( {u_{h}^{2} + u_{r}^{2} + u_{\varphi}^{2}} \right)/2}} \right\rbrack_{R} -}} \\{\left\lbrack {E_{int} + {p/\rho} + {\left( {u_{h}^{2} + u_{r}^{2} + u_{\varphi}^{2}} \right)/2}} \right\rbrack_{L}} \\{= {\left\lbrack {E_{int} + {p/\rho} + {\left( {u_{h}^{2} + u_{\varphi}^{2}} \right)/2}} \right\rbrack_{R} -}} \\{\left\lbrack {E_{int} + {p/\rho} + {\left( {u_{h}^{2} + u_{\varphi}^{2}} \right)/2}} \right\rbrack_{L}.}\end{matrix} & (46)\end{matrix}$

taking into account that (u_(r))_(R)=(u_(r))_(L), see (25). So, fullpower contributed to the flow by N blade elements located at radius requalsN·ΔP=2πrρu _(h)(H _(R) −H _(L))·Δr.  (47)

Although relationship (47) has the same form as (43), it includes fullpower P instead mechanical power P*, and uses another definition ofspecific enthalpy H.

Let's suppose that rate of heat transfer between gas and turbine bladesis much less than total energy flux through the blade wheel and outputmechanical power (if any). If this assumption is valid, then full powerwithdrawn from the flow −ΔP should be spent for producing outputmechanical power ΔW, i.e.ΔP+ΔW=0.  (48)

Substitution of expressions for ΔW (40) and ΔP (47) in equation (48)finally results in the following independent relationship representingenergy jump condition:(H _(R) −H _(L))+rΩ[(u _(φ))_(R)−(u _(φ))_(L)]=0,  (49)

where the difference (H_(R)−H_(L)) is defined by formula (46).Relationship (49) in combination with a given equation of state p=ρ(p,T) and expression for internal energy E_(int)=E_(int)(p, T) providesmathematical closure of the problem.

Numerical modeling some types of turbo machines does not allowefficiently using the cylindrical coordinate system (h, r, φ) introducedabove. Let's consider, for example a wide angle ducted turbo machineschematically shown in FIG. 5.

In such a machine flow area is bounded by two rigid surfaces, centralbody and external duct, so that in case of steady state and axiallysymmetric flow field total mass flux is the same in all cross sectionsh=const. Hence, an average axial flow velocity varies along the turbinedepending on local area of its cross section.

In the suggested disk actuator model jump condition (23) derived incylindrical coordinates provides conservation of mass flux if and onlyif radial positions of the left and right bounds of integration area Aare exactly equal. Therefore, the model generally does not preserve massconservation in case of variable cross section area because of violationof mass balance at finite thickness of disk actuator D(r). In this casethe cylindrical coordinate system (h, r, φ) is not convenient forproperly constructing disk actuator model, and a transformed coordinatesystem should be used instead.

Let's introduce an auxiliary curvilinear coordinate system (ξ, η) in alateral section φ=const, so that the inner and outer bounds of flowfield coincide with some coordinate lines η=const. If geometrical shapesof the central body and external duct are described by given functionsr=r_(min)(h) and r=r_(max)(h), respectively, as shown in FIG. 4, thenrequired coordinate transformation (ξ, η)→(h, r) can be defined as

$\begin{matrix}\left\{ \begin{matrix}{{{h\left( {\xi,\eta} \right)} = \xi},} \\{{{r\left( {\xi,\eta} \right)} = {{r_{\min}(\xi)} + {\eta \cdot \left\lbrack {{r_{\max}(\xi)} - {r_{\min}(\xi)}} \right\rbrack}}},}\end{matrix} \right. & (50)\end{matrix}$

and r(ξ, 0)=r_(min)(ξ), r(ξ, 1)=r_(max)(ξ). Navier-Stokes equations (2)describing steady state axially symmetric flow take the following formin coordinates (ξ, η):

$\begin{matrix}{{\frac{\partial\left( {E^{*} - E_{v}^{*}} \right)}{\partial\xi} + \frac{\partial\left( {G^{*} - G_{v}^{*}} \right)}{\partial\eta}} = {S^{*} - {S_{v}^{*}.}}} & (51)\end{matrix}$Transformed flux and source vectors in equations (51) are expressed interms respective vectors in equations (2) asE* _((v))=[(∂ξ/∂h)E _((v))+(∂ξ/∂r)G _((v))]·det(J),G* _((v))=[(∂η/∂h)E _((v))+(∂η/∂r)G _((v))]·det(J),S* _((v)) =S _((v))·det(J),  (52)

where J is Jacobian matrix of the coordinate transformation (ξ, η)→(h,r)

${J = {\frac{\partial\left( {r,h} \right)}{\partial\left( {\xi,\eta} \right)} = {\begin{matrix}{{\partial h}/{\partial\xi}} & {{\partial h}/{\partial\eta}} \\{{\partial r}/{\partial\xi}} & {{\partial r}/{\partial\eta}}\end{matrix}}}},{J^{- 1} = {\frac{\partial\left( {\xi,\eta} \right)}{\partial\left( {r,h} \right)} = {\begin{matrix}{{\partial\xi}/{\partial h}} & {{\partial\xi}/{\partial r}} \\{{\partial\overset{.}{\eta}}/{\partial h}} & {{\partial\eta}/{\partial r}}\end{matrix}}}}$For transformation (50) ∂h/∂ξ=1, ∂h/∂η=0,

$\begin{matrix}{{J = {\begin{matrix}1 & 0 \\r_{\xi}^{\prime} & r_{\eta}^{\prime}\end{matrix}}},{{\det\mspace{11mu}(J)} = r_{\eta}^{\prime}},{J^{- 1} = {\begin{matrix}1 & 0 \\{{- r_{\xi}^{\prime}}/r_{\eta}^{\prime}} & {1/r_{\eta}^{\prime}}\end{matrix}}}} & (53)\end{matrix}$

where brief notations for the derivatives r′_(ξ)=∂r/∂ξ and r′_(η)=∂r/∂ηare introduced. After substitution of metric coefficients (53) inexpressions (52) the latter reduce toE* _((v)) =r′ _(η) E _((v)) , G* _((v)) =G _((v)) −r′ _(ξ) E _((v)) , S*_((v)) =r′ _(η) S _((v)),  (54)

and component-wise representations of the transformed flux and sourcevectors are

$\begin{matrix}{{E^{*} = {\begin{matrix}{{rr}_{\eta}^{\prime}\rho\; u_{h}} \\{{rr}_{\eta}^{\prime}\left( {{\rho\; u_{h}^{2}} + p} \right)} \\{{rr}_{\eta}^{\prime}\rho\; u_{h}u_{r}} \\{r^{2}r_{\eta}^{\prime}\rho\; u_{h}u_{\varphi}} \\{{rr}_{\eta}^{\prime}\rho\; u_{h}H^{(\;*\;)}}\end{matrix}}},{G^{*} = {\begin{matrix}{r\;{\rho\left( {u_{r} - {u_{h}r_{\xi}^{\prime}}} \right)}} \\{r\left\lbrack {{\rho\;{u_{h}\left( {u_{r} - {u_{h}r_{\xi}^{\prime}}} \right)}} - {pr}_{\xi}^{\prime}} \right\rbrack} \\{r\left\lbrack {{\rho\;{u_{r}\left( {u_{r} - {u_{h}r_{\xi}^{\prime}}} \right)}} + p} \right\rbrack} \\{r^{2}\rho\;{u_{\varphi}\left( {u_{r} - {u_{h}r_{\xi}^{\prime}}} \right)}} \\{r\;\rho\;{H^{(\;*\;)}\left( {u_{r} - {u_{h}r_{\xi}^{\prime}}} \right)}}\end{matrix}}},{S^{*} = {\begin{matrix}0 \\{{rr}_{\eta}^{\prime}f_{h}} \\{r_{\eta}^{\prime}\left( {p + {\rho\; u_{\varphi}^{2}} + {rf}_{r}} \right)} \\{r^{2}r_{\eta}^{\prime}f_{\varphi}} \\{{rr}^{\prime}e}\end{matrix}}}} & (55)\end{matrix}$

Let's select an infinitely narrow integration area A in plane (ξ, η), sothat its lower and upper bounds cross lateral section of disk actuatorand are located on coordinate lines η and η+Δη, respectively, while itsleft (ξ=ξ_(L)(η)) and right (ξ=ξ_(R)(η)) bounds are located in the areasof free stream closely adjoining left and right sides of the disk, seeFIG. 4. Similarly to our previous considerations, flow parameters andmetric coefficients at the left and right bounds of the area A will aremarked with the indices “L” and “R”, respectively. As long as area A isdistorted in cylindrical coordinate system (h, r), radial positions ofits left and right bounds do not coincide (r_(L)≠r_(R) at Δη→0) if diskactuator of a finite thickness D(η)≠0 is used for modeling.

Vector equation (51) can be re-written in the scalar form(∂E* _(k)/∂ξ)+(∂G* _(k)/∂η)=S* _(k).  (56)

Integration of equations (56) over area A, application of the 1-stGreen's integral formula, and representation of respective integrals interms of average values (see above), result in averaged equations[E* _(k)(η)]_(R) −[E* _(k)(η)]_(L) =D[S* _(k)(η)−∂G*_(k)(η)/∂η]_(avr),  (57)

which are quite similar to equations (9). Jump conditions for flowparameters can be directly derived from equations (57) can in exactlythe same way as above:

$\quad\left\{ \begin{matrix}{{\left( {{rr}_{\eta}^{\prime}\rho\; u_{h}} \right)_{R} = \left( {{rr}_{\eta}^{\prime}\rho\; u_{h}} \right)_{L}},\mspace{500mu}(58)} \\{{{\left( u_{h} \right)_{R} - \left( u_{h} \right)_{L} + {\left( {p_{R} - p_{L}} \right){\left( {rr}_{\eta}^{\prime} \right)_{avr}/\left( {{rr}_{\eta}^{\prime}\rho\; u_{h}} \right)}}} = {D \cdot {\left( {{rr}_{\eta}^{\prime}f_{h}} \right)_{avr}/\left( {{rr}_{\eta}^{\prime}\rho\; u_{h}} \right)}}},\;(59)} \\{{\left( u_{r} \right)_{R} = \left( u_{r} \right)_{L}},\mspace{596mu}(60)} \\{{{\left( {ru}_{\varphi} \right)_{R} - \left( {ru}_{\varphi} \right)_{L}} = {D \cdot {\left( {r^{2}r_{\eta}^{\prime}f_{\varphi}} \right)_{avr}/\left( {{rr}_{\eta}^{\prime}\rho\; u_{h}} \right)}}},\mspace{284mu}(61)} \\{{H_{R}^{(\;*\;)} - H_{L}^{(\;*\;)}} = {D \cdot {\left( {{rr}_{\eta}^{\prime}e} \right)_{avr}/{\left( {{rr}_{\eta}^{\prime}\rho\; u_{h}} \right).\mspace{346mu}(62)}}}}\end{matrix} \right.$

It can be shown that approximate relationships (59) and (60) are validwith relative accuracy O(ε·δ), while (58), (61) and (62) are valid withrelative accuracy O(ε*·δ), whereε*=(u _(r) −u _(h) r′ _(ξ))/u _(h).  (63)

Note that magnitude of parameter ε* really should be small, as long asε* characterizes angular difference between directions of streamlinesand respective coordinate lines η=const. In particular, near flow boundsu_(r)/u_(h)=r′_(ξ), i.e. ε*=0 at both η=0 and η=1. On the other hand,magnitude of parameters generally is not assumed to be small, since inthis case ε=u_(r)/u_(h)=O(rr′_(ξ)/r′_(η))=O(1). Hence ε*=O(ε), andoverall relative accuracy of the model is O(ε·δ).

The Blade Element Momentum Theory can now be used for representingvolume densities of forces f_(h) and f_(φ) in terms of drag and liftcoefficients C_(D), C_(L) of blade elements. Substitution of thepreviously derived expressions (34)-(37) into (59) and (61) results inthe following representations of h- and φ-momentum jumps:[(u _(h))_(R)−(u _(h))_(L)](rr′ _(η) ρu _(h))/(rr′ _(η))_(avr)+(p _(R)−p _(L))=qσ(C _(L) sin β_(avr) −C _(D) cos β_(avr)),  (64)[(ru _(φ))_(R)−(ru _(φ))_(L)](rr′ _(η) ρu _(h))/(r ² r′ _(η))_(avr)=qσ(C _(D) sin β_(avr) +C _(L) cos β_(avr)).  (65)

Relationships (58), (60), (64), and (65) constitute a complete set ofjump conditions for equations of continuity and momentums. Jumpcondition (62) for energy equation should be considered in accordancewith the previous conclusions made in above for incompressible fluids orfor compressible gases. Total drag, torque, and mechanical powerproduced by the blade wheel are

$\begin{matrix}{{X = {{- 2}\pi{\int_{0}^{1}{\left\{ {{\left\lbrack {\left( u_{h} \right)_{R} - \left( u_{h} \right)_{L}} \right\rbrack\left( {{rr}_{\eta}^{\prime}\rho\; u_{h}} \right)} + {\left( {p_{R} - p_{L}} \right)\left( {rr}_{\eta}^{\prime} \right)_{avr}}} \right\} \cdot \ {\mathbb{d}\eta}}}}},{M = {2\pi{\int_{0}^{1}{\left\lbrack {\left( {ru}_{\varphi} \right)_{R} - \left( {ru}_{\varphi} \right)_{L}} \right\rbrack{\left( {{rr}_{\eta}^{\prime}\rho\; u_{h}} \right) \cdot \ {\mathbb{d}\eta}}}}}},{W = {\Omega \cdot M}},} & (66)\end{matrix}$

where the integrand expressions are represented by (64) and (65). Notethat the product (rr′_(η) ρu_(h)) is specifically isolated in formulas(59), (61), (62), and (64)-(66) as rr′_(η)ρu_(h)=const at η=const inaccordance with (58). It can be seen that relationships (64), (65), and(66) are quite similar to the previously derived (38), (39), and (41),respectively.

Jump conditions (58)-(62) can be used instead of conditions (23)-(27)not only for keeping mass balance when modeling wide angle ducted turbomachines, but also in context of other similar problems, if for example,cylindrical coordinate system (h, r, φ) is not convenient for someproblem-specific reasons, or if accuracy of (23)-(27) is not sufficientbecause of a strong radial flow.

Basic conditions of applicability of the suggested model can besummarized as follows. Flow field in the area of blade wheel locationcan be considered as steady-state and axially symmetric, i.e. impact ofunsteady and azimuth variation effects on momentum and energy balancesis negligible. Scale factor (ε·δ) is small in the regions of blade wheelswept area providing major contribution to mass, momentum, and energybalances. If a flow of compressible gas is considered, then rate ofdissipative heat transfer between the gas and turbine blades is muchless than total energy flux through the blade wheel and outputmechanical power (if any).

The following references were mentioned above. 1. A. Betz. Das Maximumdes theoretisch möglichen Ausnützung des Windes durch Windmotoren,Zeitschrifft für das gesamte Turbinewesen, Volume 26, p. 307, 1920. 2.H. Glauert. Windmills and Fans. In W. F. Durand (ed). AerodynamicTheory. Dover Publications Inc., New York 1963. 3. J. Laursen, P.Enevoldsen, and S. Hjort. 3D CFD Quantification of the Performance of aMulti-Megawatt Wind Turbine. The Science of Making Torque from Wind,Journal of Physics: Conference Series 75. IOP Publishing Ltd., 2007. 4.D. Hartwanger, A. Horvat. 3D Modeling of Wind Turbine Using CFD. NAFEMSConference 2008, United Kingdom, June 2008. 5. S. S. A. Ivanell.Numerical Computations of Wind Turbine Wakes. Technical Reports fromRoyal Institute of Technology Linnre Flow Centre, Department ofMechanics, Stockholm, Sweden, January 2009. 6. L. Battisti, G.Soraperra, R. Fedrizzi, and L. Zanne. Inverse Design-Momentum, a Methodfor the Preliminary Design of Horizontal Axis Wind Turbines. The Scienceof Making Torque from Wind, Journal of Physics: Conference Series 75.IOP Publishing Ltd., 2007. 7. R. Mikkelsen. Actuator Disk MethodsApplied to Wind Turbines. Dissertation submitted to Technical Universityof Denmark, Fluid Mechanics, Department of Mechanical Engineering, 2003.

Embodiments of the subject matter and the functional operationsdescribed in this specification can be implemented in digital electroniccircuitry, in tangibly-embodied computer software or firmware, incomputer hardware, including the structures disclosed in thisspecification and their structural equivalents, or in combinations ofone or more of them. Embodiments of the subject matter described in thisspecification can be implemented as one or more computer programs, i.e.,one or more modules of computer program instructions encoded on atangible non-transitory program carrier for execution by, or to controlthe operation of, data processing apparatus. Alternatively or inaddition, the program instructions can be encoded on anartificially-generated propagated signal, e.g., a machine-generatedelectrical, optical, or electromagnetic signal, that is generated toencode information for transmission to suitable receiver apparatus forexecution by a data processing apparatus. The computer storage mediumcan be a machine-readable storage device, a machine-readable storagesubstrate, a random or serial access memory device, or a combination ofone or more of them. The computer storage medium is not, however, apropagated signal.

The term “data processing apparatus” encompasses all kinds of apparatus,devices, and machines for processing data, including by way of example aprogrammable processor, a computer, or multiple processors or computers.The apparatus can include special purpose logic circuitry, e.g., an FPGA(field programmable gate array) or an ASIC (application-specificintegrated circuit). The apparatus can also include, in addition tohardware, code that creates an execution environment for the computerprogram in question, e.g., code that constitutes processor firmware, aprotocol stack, a database management system, an operating system, or acombination of one or more of them.

A computer program (which may also be referred to or described as aprogram, software, a software application, a module, a software module,a script, or code) can be written in any form of programming language,including compiled or interpreted languages, or declarative orprocedural languages, and it can be deployed in any form, including as astand-alone program or as a module, component, subroutine, or other unitsuitable for use in a computing environment. A computer program may, butneed not, correspond to a file in a file system. A program can be storedin a portion of a file that holds other programs or data, e.g., one ormore scripts stored in a markup language document, in a single filededicated to the program in question, or in multiple coordinated files,e.g., files that store one or more modules, sub-programs, or portions ofcode. A computer program can be deployed to be executed on one computeror on multiple computers that are located at one site or distributedacross multiple sites and interconnected by a communication network.

As used in this specification, an “engine,” or “software engine,” refersto a software implemented input/output system that provides an outputthat is different from the input. An engine can be an encoded block offunctionality, such as a library, a platform, a software development kit(“SDK”), or an object. Each engine can be implemented on any appropriatetype of computing device, e.g., servers, mobile phones, tabletcomputers, notebook computers, music players, e-book readers, laptop ordesktop computers, PDAs, smart phones, or other stationary or portabledevices, that includes one or more processors and computer readablemedia. Additionally, two or more of the engines may be implemented onthe same computing device, or on different computing devices.

The processes and logic flows described in this specification can beperformed by one or more programmable computers executing one or morecomputer programs to perform functions by operating on input data andgenerating output. The processes and logic flows can also be performedby, and apparatus can also be implemented as, special purpose logiccircuitry, e.g., an FPGA (field programmable gate array) or an ASIC(application-specific integrated circuit).

Computers suitable for the execution of a computer program include, byway of example, can be based on general or special purposemicroprocessors or both, or any other kind of central processing unit.Generally, a central processing unit will receive instructions and datafrom a read-only memory or a random access memory or both. The essentialelements of a computer are a central processing unit for performing orexecuting instructions and one or more memory devices for storinginstructions and data. Generally, a computer will also include, or beoperatively coupled to receive data from or transfer data to, or both,one or more mass storage devices for storing data, e.g., magnetic,magneto-optical disks, or optical disks. However, a computer need nothave such devices. Moreover, a computer can be embedded in anotherdevice, e.g., a mobile telephone, a personal digital assistant (PDA), amobile audio or video player, a game console, a Global PositioningSystem (GPS) receiver, or a portable storage device, e.g., a universalserial bus (USB) flash drive, to name just a few.

Computer-readable media suitable for storing computer programinstructions and data include all forms of non-volatile memory, mediaand memory devices, including by way of example semiconductor memorydevices, e.g., EPROM, EEPROM, and flash memory devices; magnetic disks,e.g., internal hard disks or removable disks; magneto-optical disks; andCD-ROM and DVD-ROM disks. The processor and the memory can besupplemented by, or incorporated in, special purpose logic circuitry.

To provide for interaction with a user, embodiments of the subjectmatter described in this specification can be implemented on a computerhaving a display device, e.g., a CRT (cathode ray tube) or LCD (liquidcrystal display) monitor, for displaying information to the user and akeyboard and a pointing device, e.g., a mouse or a trackball, by whichthe user can provide input to the computer. Other kinds of devices canbe used to provide for interaction with a user as well; for example,feedback provided to the user can be any form of sensory feedback, e.g.,visual feedback, auditory feedback, or tactile feedback; and input fromthe user can be received in any form, including acoustic, speech, ortactile input. In addition, a computer can interact with a user bysending documents to and receiving documents from a device that is usedby the user; for example, by sending web pages to a web browser on auser's client device in response to requests received from the webbrowser.

Embodiments of the subject matter described in this specification can beimplemented in a computing system that includes a back-end component,e.g., as a data server, or that includes a middleware component, e.g.,an application server, or that includes a front-end component, e.g., aclient computer having a graphical user interface or a Web browserthrough which a user can interact with an implementation of the subjectmatter described in this specification, or any combination of one ormore such back-end, middleware, or front-end components. The componentsof the system can be interconnected by any form or medium of digitaldata communication, e.g., a communication network. Examples ofcommunication networks include a local area network (“LAN”) and a widearea network (“WAN”), e.g., the Internet.

The computing system can include clients and servers. A client andserver are generally remote from each other and typically interactthrough a communication network. The relationship of client and serverarises by virtue of computer programs running on the respectivecomputers and having a client-server relationship to each other.

While this specification contains many specific implementation details,these should not be construed as limitations on the scope of anyinvention or of what may be claimed, but rather as descriptions offeatures that may be specific to particular embodiments of particularinventions. Certain features that are described in this specification inthe context of separate embodiments can also be implemented incombination in a single embodiment. Conversely, various features thatare described in the context of a single embodiment can also beimplemented in multiple embodiments separately or in any suitablesubcombination. Moreover, although features may be described above asacting in certain combinations and even initially claimed as such, oneor more features from a claimed combination can in some cases be excisedfrom the combination, and the claimed combination may be directed to asubcombination or variation of a subcombination.

Similarly, while operations are depicted in the drawings in a particularorder, this should not be understood as requiring that such operationsbe performed in the particular order shown or in sequential order, orthat all illustrated operations be performed, to achieve desirableresults. In certain circumstances, multitasking and parallel processingmay be advantageous. Moreover, the separation of various system modulesand components in the embodiments described above should not beunderstood as requiring such separation in all embodiments, and itshould be understood that the described program components and systemscan generally be integrated together in a single software product orpackaged into multiple software products.

Particular embodiments of the subject matter have been described. Otherembodiments are within the scope of the following claims. For example,the actions recited in the claims can be performed in a different orderand still achieve desirable results. As one example, the processesdepicted in the accompanying figures do not necessarily require theparticular order shown, or sequential order, to achieve desirableresults. In certain implementations, multitasking and parallelprocessing may be advantageous.

What is claimed is:
 1. A computer-implemented method comprising:obtaining a plurality of parameters of a design of a blade wheel of anaxial flow turbo machine, including lift coefficients and dragcoefficients for the blade wheel for each of a plurality of angles ofattack α_(i) and for each of a plurality of cross sections r_(j) of theblade wheel; from a minimum radius r_(min) of the blade wheel to amaximum radius r_(max) of the blade wheel; obtaining mass jumpconditions, momentum jump conditions, and energy jump conditions of atwo-dimensional blade element model of the blade wheel, wherein the massjump conditions, momentum jump conditions, and energy jump conditionseach specify a relationship between upstream flow field parameters anddownstream flow field parameters of the blade wheel; computing upstreamflow field parameters and downstream flow field parameters that satisfythe mass jump conditions, momentum jump conditions, and energy jumpconditions of the two-dimensional blade element model; computingmomentum jump values across the swept area of the blade wheel using thecomputed upstream flow field parameters and the computed downstream flowfield parameters that satisfy the mass jump conditions, momentum jumpconditions, and energy jump conditions of the two-dimensional bladeelement model; computing one or more turbine parameters for the axialflow turbo machine including computing an approximation of total torquefor the design of the blade wheel by aggregating the momentum jumpvalues across the swept area of the blade wheel; and outputting thedesign of the blade wheel based on the computed one or more turbineparameters.
 2. The method of claim 1, wherein computing the momentumjump values from the computed upstream flow field parameters and thecomputed downstream flow field parameters comprises: computing averageddirections β_(i) of velocity components of the upstream flow fieldparameters and the downstream flow field parameters; and computing eachmomentum jump values m_(i) according to: m_(i)=q·σ·[C_(D) sinβ_(i)+C_(L) cos β_(i)], where q is an average dynamic pressure, C_(D)are drag coefficients of the blade wheel, C_(L) are lift coefficients ofthe blade wheel, and wherein σ is a solidity factor of the blade wheel.3. The method of claim 1, wherein computing the upstream flow fieldparameters and the downstream flow field parameters comprises computingupstream density ρ_(L), upstream pressure p_(L), upstream axial velocity(u_(h))_(L), upstream azimuth velocity (u_(φ))_(L), downstream densityρ_(R), downstream pressure p_(R), downstream axial velocity (u_(h))_(R),and downstream azimuth velocity (u_(φ))_(R), and further comprising:computing an approximation of total drag X produced by the blade wheelaccording to:X = −2π∫_(r_(min))^(r_(max))r{ρ u_(h)[(u_(h))_(R) − (u_(h))_(L)] + (p_(R) − p_(L))} ⋅ 𝕕r.4. The method of claim 1, wherein computing the upstream flow fieldparameters and the downstream flow field parameters comprises computingupstream density ρ_(L), upstream pressure p_(L), upstream axial velocity(u_(h))_(L), upstream azimuth velocity (u_(φ))_(L), downstream densityρ_(R), downstream pressure p_(R), downstream axial velocity (u_(h))_(R),and downstream azimuth velocity (u_(φ))_(R), and further comprising:computing the approximation of the total torque M according to:M = 2π∫_(r_(min))^(r_(max))r²ρ u_(h)[(u_(φ))_(R) − (u_(φ))_(L)] ⋅ 𝕕r. 5.The method of claim 4, further comprising computing an approximation oftotal mechanical power W according to:W=Ω·M, wherein Ω is the rotational speed of the blade wheel.
 6. Themethod of claim 4, wherein computing the upstream flow field parametersand the downstream flow field parameters comprises computing upstreamtotal enthalpy of the blade wheel H*_(L) and downstream total enthalpyof the blade wheel H*_(R) and further comprising: computing anapproximation of total mechanical power P* withdrawn from flow by theblade wheels in an incompressible fluid according to:P^(*) = 2π∫_(r_(min))^(r_(max))r ρ u_(h)(H_(R)^(*) − H_(L)^(*)) ⋅ 𝕕r. 7.The method of claim 1, wherein the turbo machine is a wide angle ductedturbo machine, and further comprising numerically simulating theupstream flow field parameters and the downstream flow field parametersusing representations of radial distributions of mass, momentum, andenergy jumps across the swept area of the blade wheel derived withboundary-fitted coordinate transformation h(ξ, η) for coordinates ξ andη as h(ξ, η)=ξ, r(ξ, η)=r_(min)(ξ)+η·[r_(max)(ξ)−r_(min)(ξ)].
 8. Asystem comprising: one or more computers and one or more storage devicesstoring instructions that are operable, when executed by the one or morecomputers, to cause the one or more computers to perform operationscomprising: obtaining a plurality of parameters of a design of a bladewheel of an axial flow turbo machine, including lift coefficients anddrag coefficients for the blade wheel for each of a plurality of anglesof attack α_(i) and for each of a plurality of cross sections r_(j) ofthe blade wheel from a minimum radius r_(min) of the blade wheel to amaximum radius r_(max) of the blade wheel; obtaining mass jumpconditions, momentum jump conditions, and energy jump conditions of atwo-dimensional blade element model of the blade wheel, wherein the massjump conditions, momentum jump conditions, and energy jump conditionseach specify a relationship between upstream flow field parameters anddownstream flow field parameters of the blade wheel; computing upstreamflow field parameters and downstream flow field parameters that satisfythe mass jump conditions, momentum jump conditions, and energy jumpconditions of the two-dimensional blade element model; computingmomentum jump values across the swept area of the blade wheel using thecomputed upstream flow field parameters and the computed downstream flowfield parameters that satisfy the mass jump conditions, momentum jumpconditions, and energy jump conditions of the two-dimensional bladeelement model; and computing one or more turbine parameters for theaxial flow turbo machine including computing an approximation of totaltorque for the design of the blade wheel by aggregating the momentumjump values across the swept area of the blade wheel; and outputting thedesign of the blade wheel based on the computed one or more turbineparameters.
 9. The system of claim 8, wherein computing the momentumjump values from the computed upstream flow field parameters and thecomputed downstream flow field parameters comprises: computing averageddirections β_(i) of velocity components of the upstream flow fieldparameters and the downstream flow field parameters; and computing eachmomentum jump values m_(i) according to: m_(i)=q·σ·[C_(D) sinβ_(i)+C_(L) cos β_(i)], where q is an average dynamic pressure, C_(D)are drag coefficients of the blade wheel, C_(L) are lift coefficients ofthe blade wheel, and wherein σ is a solidity factor of the blade wheel.10. The system of claim 8, wherein computing the upstream flow fieldparameters and the downstream flow field parameters comprises computingupstream density ρ_(L), upstream pressure p_(L), upstream axial velocity(u_(h))_(L), upstream azimuth velocity (u_(φ))_(L), downstream densityρ_(R), downstream pressure p_(R), downstream axial velocity (u_(h))_(R),and downstream azimuth velocity (u_(φ))_(R), and wherein the operationsfurther comprise: computing an approximation of total drag X produced bythe blade wheel according to:X = −2π∫_(r_(min))^(r_(max))r{ρ u_(h)[(u_(h))_(R) − (u_(h))_(L)] + (p_(R) − p_(L))} ⋅ 𝕕r.11. The system of claim 8, wherein computing the upstream flow fieldparameters and the downstream flow field parameters comprises computingupstream density ρ_(L), upstream pressure p_(L), upstream axial velocity(u_(h))_(L), upstream azimuth velocity (u_(φ))_(L), downstream densityρ_(R), downstream pressure p_(R), downstream axial velocity (u_(h))_(R),and downstream azimuth velocity (u_(φ))_(R), and wherein the operationsfurther comprise: computing the approximation of the total torque Maccording to:M = 2π∫_(r_(min))^(r_(max))r²ρ u_(h)[(u_(φ))_(R) − (u_(φ))_(L)] ⋅ 𝕕r.12. The system of claim 11, wherein the operations further comprisecomputing an approximation of total mechanical power W according to:W=Ω·M, wherein Ω is the rotational speed of the blade wheel.
 13. Thesystem of claim 11, wherein computing the upstream flow field parametersand the downstream flow field parameters comprises computing upstreamtotal enthalpy of the blade wheel H*_(L) and downstream total enthalpyof the blade wheel H*_(R) and wherein the operations further comprise:computing an approximation of total mechanical power P* withdrawn fromflow by the blade wheels in an incompressible fluid according to:P^(*) = 2π∫_(r_(min))^(r_(max))r ρ u_(h)(H_(R)^(*) − H_(L)^(*)) ⋅ 𝕕r.14. The system of claim 8, wherein the turbo machine is a wide angleducted turbo machine and wherein the operations further comprise:numerically simulating the upstream flow field parameters and thedownstream flow field parameters using representations of radialdistributions of mass, momentum, and energy jumps across the swept areaof the blade wheel derived with boundary-fitted coordinatetransformation h(ξ, η) for coordinates ξ and η as h(ξ, η)=ξ, r(ξ,η)=r_(min)(ξ)+η·[r_(max)(ξ)−r_(min)(ξ)].
 15. A computer program product,encoded on one or more non-transitory computer storage media, comprisinginstructions that when executed by one or more computers cause the oneor more computers to perform operations comprising: obtaining aplurality of parameters of a design of a blade wheel of an axial flowturbo machine, including lift coefficients and drag coefficients for theblade wheel for each of a plurality of angles of attack α_(i) and foreach of a plurality of cross sections r_(j) of the blade wheel from aminimum radius r_(min) of the blade wheel to a maximum radius r_(max) ofthe blade wheel; obtaining mass jump conditions, momentum jumpconditions, and energy jump conditions of a two-dimensional bladeelement model of the blade wheel, wherein the mass jump conditions,momentum jump conditions, and energy jump conditions each specify arelationship between upstream flow field parameters and downstream flowfield parameters of the blade wheel; computing upstream flow fieldparameters and downstream flow field parameters that satisfy the massjump conditions, momentum jump conditions, and energy jump conditions ofthe two-dimensional blade element model; computing momentum jump valuesacross the swept area of the blade wheel using the computed upstreamflow field parameters and the computed downstream flow field parametersthat satisfy the mass jump conditions, momentum jump conditions, andenergy jump conditions of the two-dimensional blade element model; andcomputing one or more turbine parameters for the axial flow turbomachine including computing an approximation of total torque for thedesign of the blade wheel by aggregating the momentum jump values acrossthe swept area of the blade wheel; and outputting the design of theblade wheel based on the computed one or more turbine parameters. 16.The computer program product of claim 15, wherein computing the momentumjump values from the computed upstream flow field parameters and thecomputed downstream flow field parameters comprises: computing averageddirections β_(i), of velocity components of the upstream flow fieldparameters and the downstream flow field parameters; and computing eachmomentum jump values m_(i) according to: m_(i)=q·σ·[C_(D) sinβ_(i)+C_(L) cos β_(i)], where q is an average dynamic pressure C_(D) aredrag coefficients of the blade wheel, C_(L) are lift coefficients of theblade wheel, and wherein σ is a solidity factor of the blade wheel. 17.The computer program product of claim 15, wherein computing the upstreamflow field parameters and the downstream flow field parameters comprisescomputing upstream density ρ_(L), upstream pressure p_(L), upstreamaxial velocity (u_(h))_(L), upstream azimuth velocity (u_(φ))_(L),downstream density ρ_(R), downstream pressure p_(R), downstream axialvelocity (u_(h))_(R), and downstream azimuth velocity (u_(φ))_(R), andwherein the operations further comprise: computing an approximation oftotal drag X produced by the blade wheel according to:X = −2π∫_(r_(min))^(r_(max))r{ρ u_(h)[(u_(h))_(R) − (u_(h))_(L)] + (p_(R) − p_(L))} ⋅ 𝕕r.18. The computer program product of claim 15, wherein computing theupstream flow field parameters and the downstream flow field parameterscomprises computing upstream density ρ_(L), upstream pressure p_(L),upstream axial velocity (u_(h))_(L), upstream azimuth velocity(u_(φ))_(L), downstream density ρ_(R), downstream pressure p_(R),downstream axial velocity (u_(h))_(R), and downstream azimuth velocity(u_(φ))_(R), and wherein the operations further comprise: computing theapproximation of the total torque M according to:M = 2π∫_(r_(min))^(r_(max))r²ρ u_(h)[(u_(φ))_(R) − (u_(φ))_(L)] ⋅ 𝕕r.19. The computer program product of claim 18, wherein the operationsfurther comprise computing an approximation of total mechanical power Waccording to:W=Ω·M, wherein Ω is the rotational speed of the blade wheel.
 20. Thecomputer program product of claim 18, wherein computing the upstreamflow field parameters and the downstream flow field parameters comprisescomputing upstream total enthalpy of the blade wheel H*_(L) anddownstream total enthalpy of the blade wheel H*_(R) and wherein theoperations further comprise: computing total mechanical power P*withdrawn from flow by the blade wheels in an incompressible fluidaccording to:P^(*) = 2π∫_(r_(min))^(r_(max))r ρ u_(h)(H_(R)^(*) − H_(L)^(*)) ⋅ 𝕕r.21. The computer program product of claim 15, wherein the turbo machineis a wide angle ducted turbo machine and wherein the operations furthercomprise: numerically simulating the upstream flow field parameters andthe downstream flow field parameters using representations of radialdistributions of mass, momentum, and energy jumps across the swept areaof the blade wheel derived with boundary-fitted coordinatetransformation h(ξ, η) for coordinates ξ and η as h(ξ, η)=ξ, r(ξ,η)=r_(min)(ξ)+η·[r_(max)(ξ)−r_(min)(ξ)].